##Data for Figure 1 and Figure 2
library(foreign)
## YOU MUST SET YOUR WORKING DIRECTORY HERE USING SETWD() OR HERE() 
COVID=read.dta("COVIDdataforR.dta")

##Figure 1
COVID.R=subset(COVID,COVID$republican==1)
COVID.D=subset(COVID,COVID$republican==0)

WORRY.R=density(COVID.R$worry,na.rm=T)
WORRY.D=density(COVID.D$worry,na.rm=T)
WORRY.ALL=density(COVID$worry,na.rm=T)

BEHAVIOR.R=density(COVID.R$covidbehavior,na.rm=T)
BEHAVIOR.D=density(COVID.D$covidbehavior,na.rm=T)
BEHAVIOR.ALL=density(COVID$covidbehavior,na.rm=T)

POLICY.R=density(COVID.R$covidpolicy,na.rm=T)
POLICY.D=density(COVID.D$covidpolicy,na.rm=T)
POLICY.ALL=density(COVID$covidpolicy,na.rm=T)

pdf(file="figure1.pdf", 
    height = 8.5, 
    width = 11) 
par(oma=c(7,7,0,0),mfrow=c(1,3),mar=c(1,2,1,2))
plot(WORRY.R,main=" ",ylab=" ",xlab=" ",type="n",cex.axis=1.7)
lines(WORRY.R,col="red",lty=1,lwd=3)
lines(WORRY.D,col="blue",lty=1,lwd=3)

mtext("Density",side=2,cex=1,line=5)
mtext("Worry about COVID",side=1,cex=1,line=5)

plot(BEHAVIOR.D,main=" ",ylab=" ",xlab=" ",type="n",cex.axis=1.7)
lines(BEHAVIOR.R,col="red",lty=1,lwd=3)
lines(BEHAVIOR.D,col="blue",lty=1,lwd=3)

mtext("Changes in Behavior",side=1,cex=1,line=5)

plot(POLICY.D,main=" ",ylab=" ",xlab=" ",type="n",cex.axis=1.7)
lines(POLICY.R,col="red",lty=1,lwd=3)
lines(POLICY.D,col="blue",lty=1,lwd=3)

legend(-0.1,3,legend=c("Democrats","Republicans"),lwd=3,lty=1,col=c("blue","red"),cex=2,bty="n",horiz=F)
mtext("Support for COVID Policies",side=1,cex=1,line=5)
dev.off() 

##Figure 2

pdf(file="figure2.pdf",
    height = 8.5, 
    width = 11)
COVID.WORRY=subset(COVID,COVID$worry>=0)
COVID.BEHAV=subset(COVID,COVID$covidbehavior>=0)
COVID.POLICY=subset(COVID,COVID$covidpolicy>=0)
COVID.WORRY=subset(COVID.WORRY,COVID.WORRY$day3_avg_cases>=0)
COVID.BEHAV=subset(COVID.BEHAV,COVID.BEHAV$day3_avg_cases>=0)
COVID.POLICY=subset(COVID.POLICY,COVID.POLICY$day3_avg_cases>=0)

WORRY.LO=lowess(COVID.WORRY$day3_avg_cases,COVID.WORRY$worry,f=.5)
BEHAVIOR.LO=lowess(COVID.BEHAV$day3_avg_cases,COVID.BEHAV$covidbehavior,f=.5)
POLICY.LO=lowess(COVID.POLICY$day3_avg_cases,COVID.POLICY$covidpolicy,f=.5)


par(oma=c(7,6,0,0),mfrow=c(1,3),mar=c(1,5,1,5))
plot(COVID$day3_avg_cases,COVID$worry,type="n",xlim=c(0,100000),axes=F,xlab=" ",ylab=" ")
points(jitter(COVID$day3_avg_cases,factor=2),jitter(COVID$worry,factor=2),col="grey55",bg="grey55",pch=21)
lines(WORRY.LO,lty=1,lwd=3,col="black")
mtext("Worry about COVID",side=2,cex=1.5,line=5)
axis(1,at = seq(0,100000,25000), label = seq(0,100000,25000), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75,las=2)


plot(COVID$day3_avg_cases,COVID$covidbehavior,type="n",xlim=c(0,100000),axes=F,xlab=" ",ylab=" ")
points(jitter(COVID$day3_avg_cases,factor=2),jitter(COVID$covidbehavior,factor=2),col="grey55",bg="grey55",pch=21)
lines(BEHAVIOR.LO,lty=1,lwd=3,col="black")
mtext("Changes in Behaviour",side=2,cex=1.5,line=5)
axis(1,at = seq(0,100000,25000), label = seq(0,100000,25000), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(0,14,2), label = seq(0,14,2), mgp = c(.8,2,1), cex.axis=1.75,las=2)


plot(COVID$day3_avg_cases,COVID$covidpolicy,type="n",xlim=c(0,100000),axes=F,xlab=" ",ylab=" ")
points(jitter(COVID$day3_avg_cases,factor=2),jitter(COVID$covidpolicy,factor=2),col="grey55",bg="grey55",pch=21)
lines(POLICY.LO,lty=1,lwd=3,col="black")
mtext("Support for COVID Policies",side=2,cex=1.5,line=5)
axis(1,at = seq(0,100000,25000), label = seq(0,100000,25000), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("Number of Cases in the County",side=1,cex=1.5,line=3,outer=T)
dev.off() 

##Figure 3 

pdf(file="figure3.pdf",
    height = 11, 
    width = 8.5)
## Panel A 
 worrypol.me=c(0.01888,-0.13816)
 worrypol.se=c(0.02793,0.04774)

 behavpol.me=c(1.45357,-0.46282)
 behavpol.se=c(0.46065,0.60567)

policypol.me=c(0.06625,-0.18150)
policypol.se=c(0.03022,0.07207)

ruler=c(1,2)

par(mfrow=c(2,3),oma=c(4,3.5,2,0.5),mar=c(4,3.5,2,3))
plot(ruler,worrypol.me,ylim=c(-.3,.3),xlim=c(.5,2.5),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,worrypol.me-1.96*worrypol.se, 
         ruler,worrypol.me+1.96*worrypol.se, lwd =  1)
segments(ruler,worrypol.me-1.64*worrypol.se, 
         ruler,worrypol.me+1.64*worrypol.se, lwd =  3)
points(ruler, worrypol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,2,by=1), label = c("Democrat","Republican"), las = 1, tick = T, cex.axis =1.2)
axis(2, at = seq(-.3,.3,by=.1), label = c(-.3,-.2,-.1,0,.1,.2,.3), las = 1, tick = T, cex.axis =1)
mtext("Worry about COVID",line=1,side=3,cex=1)


plot(ruler,behavpol.me,ylim=c(-3,3),xlim=c(.5,2.5),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,behavpol.me-1.96*behavpol.se, 
         ruler,behavpol.me+1.96*behavpol.se, lwd =  1)
segments(ruler,behavpol.me-1.64*behavpol.se, 
         ruler,behavpol.me+1.64*behavpol.se, lwd =  3)
points(ruler, behavpol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,2,by=1), label = c("Democrat","Republican"), las = 1, tick = T, cex.axis =1.2)

axis(2, at = seq(-3,3,by=1), label = seq(-3,3,by=1), las = 1, tick = T, cex.axis =1)
mtext("Changes in Behavior",line=1,side=3,cex=1)

plot(ruler,policypol.me,ylim=c(-.4,.4),xlim=c(.5,2.5),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,policypol.me-1.96*policypol.se, 
         ruler,policypol.me+1.96*policypol.se, lwd =  1)
segments(ruler,policypol.me-1.64*policypol.se, 
         ruler,policypol.me+1.64*policypol.se, lwd =  3)
points(ruler, policypol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,2,by=1), label = c("Democrat","Republican"), las = 1, tick = T, cex.axis =1.2)
axis(2, at = seq(-.4,.4,by=.2), label = seq(-.4,.4,by=.2), las = 1, tick = T, cex.axis =1)

#mtext("Marginal Effect of Partisan Animosity",line=1,side=2,cex=1, outer=TRUE) 
mtext("     Marginal Effect of Republican                                                           Marginal Effect of Partisan Animosity ", line = 1, side = 2, cex = 1, outer = T)
mtext("Support for COVID Policies",line=1,side=3,cex=1)
#dev.off() 

##Figure 3b
covidsimp=read.table("COVID Simple.txt",header=T)

outgroup=seq(0,1,by=.05)

#pdf(file="figure3b.pdf",
#    height = 8.5,
#    width = 11)
#par(mfrow=c(1,3),oma=c(10,10,4,1),mar=c(1.5,3.5,1.5,3.5))
plot(outgroup,covidsimp$worry, ylim=c(-.2,.2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$worry.ub95,rev(covidsimp$worry.lb95)),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$worry.ub90,rev(covidsimp$worry.lb90)),col="grey80",border="NA")
lines(outgroup,covidsimp$worry,lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.2)
axis(2,at = seq(-.2,.2,.1), label = seq(-.2,.2,.1), mgp = c(.8,2,1), cex.axis=1.2,las=2)
#mtext("Worry about COVID",line=-1,side=3,cex=1)

plot(outgroup,covidsimp$behav, ylim=c(-2.5,2.5), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$behav.ub95,rev(covidsimp$behav.lb95)),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$behav.ub90,rev(covidsimp$behav.lb90)),col="grey80",border="NA")
lines(outgroup,covidsimp$behav,lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.2)
axis(2,at = seq(-2.5,2.5,.5), seq(-2.5,2.5,.5), mgp = c(.8,2,1), cex.axis=1.2,las=2)
#mtext("Changes in Behavior",line=-1,side=3,cex=1)

plot(outgroup,covidsimp$policy, ylim=c(-.3,.3), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$policy.ub95,rev(covidsimp$policy.lb95)),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covidsimp$policy.ub90,rev(covidsimp$policy.lb90)),col="grey80",border="NA")
lines(outgroup,covidsimp$policy,lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.2)
axis(2,at = seq(-.3,.3,.1), label = c(-.3,-.2,-.1,0,.1,.2,.3), mgp = c(.8,2,1), cex.axis=1.2,las=2)
#mtext("Support for COVID Policies",line=-1,side=3,cex=1)
#mtext("Marginal Effect of Republican",side=2,line=2,cex=1,outer=T)
#mtext("on Dependent Variable",side=2,line=1.75,cex=1,outer=T)
mtext("Partisan Animosity",side=1,line=1,cex=1,outer=T)
dev.off() 


##Values for Figures 4 to 6
covid=read.table("COVID Triple.txt",header=T)

outgroup=seq(0,1,by=.05)

##Low 25th precentile 53 Cases (e.g., Scott County, Iowa on April 7)
##High 75th precentile 1,177 Cases (e.g., Erie County, New York on April 8)

##Figure 4 WORRY
pdf(file="figure4.pdf",
    height = 8.5,
    width = 11)
par(mfrow=c(1,2),oma=c(7,8,3,1),mar=c(1,2,1,1))
plot(outgroup,covid$worry[1:21], ylim=c(-.2,.2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$worry.ub95[1:21],rev(covid$worry.lb95[1:21])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$worry.ub90[1:21],rev(covid$worry.lb90[1:21])),col="grey80",border="NA")
lines(outgroup,covid$worry[1:21],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-.2,.2,.1), label = seq(-.2,.2,.1), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("25th precentile 53 Cases",side=3,line=1.5,cex=1.25,outer=F)
mtext("(e.g., Scott County, Iowa on April 7)",side=3,line=-0.3,cex=1,outer=F)

plot(outgroup,covid$worry[22:42], ylim=c(-.2,.2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$worry.ub95[22:42],rev(covid$worry.lb95[22:42])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$worry.ub90[22:42],rev(covid$worry.lb90[22:42])),col="grey80",border="NA")
lines(outgroup,covid$worry[22:42],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-.2,.2,.1), label = rep(" ",times=5), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("75th precentile 1,177 Cases",side=3,line=1.4,cex=1.25,outer=F)
mtext("(e.g., Erie County, New York on April 8)",side=3,line=-0.3,cex=1,outer=F)

mtext("Partisan Animosity",side=1,line=4,cex=1.25,outer=T)
mtext("Marginal Effect of ",side=2,line=5,cex=1.25,outer=T)
mtext("Republican on Worry",side=2,line=3,cex=1.25,outer=T)
dev.off() 

##Figure 5 BEHAVIORS 
pdf(file = "figure5.pdf",
    height = 8.5,
    width = 11)
par(mfrow=c(1,2),oma=c(7,8,3,1),mar=c(1,2,1,1))
plot(outgroup,covid$behav[1:21], ylim=c(-2.5,2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$behav.ub95[1:21],rev(covid$behav.lb95[1:21])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$behav.ub90[1:21],rev(covid$behav.lb90[1:21])),col="grey80",border="NA")
lines(outgroup,covid$behav[1:21],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-2.5,2,.5), seq(-2.5,2,.5), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("25th precentile 53 Cases",side=3,line=1.5,cex=1.25,outer=F)
mtext("(e.g., Scott County, Iowa on April 7)",side=3,line=-0.3,cex=1,outer=F)

plot(outgroup,covid$behav[22:42], ylim=c(-2.5,2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$behav.ub95[22:42],rev(covid$behav.lb95[22:42])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$behav.ub90[22:42],rev(covid$behav.lb90[22:42])),col="grey80",border="NA")
lines(outgroup,covid$behav[22:42],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-2.5,2,.5), label = rep(" ",times=10), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("75th precentile 1,177 Cases",side=3,line=1.4,cex=1.25,outer=F)
mtext("(e.g., Erie County, New York on April 8)",side=3,line=-0.3,cex=1,outer=F)

mtext("Partisan Animosity",side=1,line=4,cex=1.25,outer=T)
mtext("Marginal Effect of ",side=2,line=5,cex=1.25,outer=T)
mtext("Republican on Behaviors",side=2,line=3,cex=1.25,outer=T)
dev.off() 

##Figure 6 POLICY 
pdf(file = "figure6.pdf",
    height = 8.5, 
    width = 11)
par(mfrow=c(1,2),oma=c(7,8,3,1),mar=c(1,2,1,1))
plot(outgroup,covid$policy[1:21], ylim=c(-.2,.2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$policy.ub95[1:21],rev(covid$policy.lb95[1:21])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$policy.ub90[1:21],rev(covid$policy.lb90[1:21])),col="grey80",border="NA")
lines(outgroup,covid$policy[1:21],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-.2,.2,.1), seq(-.2,.2,.1), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("25th precentile 53 Cases",side=3,line=1.5,cex=1.25,outer=F)
mtext("(e.g., Scott County, Iowa on April 7)",side=3,line=-0.3,cex=1,outer=F)

plot(outgroup,covid$policy[22:42], ylim=c(-.2,.2), ylab=" ",xlab=" ",axes=F,type="n")
polygon(c(outgroup,rev(outgroup)),c(covid$policy.ub95[22:42],rev(covid$policy.lb95[22:42])),col="grey50",border="NA")
polygon(c(outgroup,rev(outgroup)),c(covid$policy.ub90[22:42],rev(covid$policy.lb90[22:42])),col="grey80",border="NA")
lines(outgroup,covid$policy[22:42],lwd=2,col="black")

abline(0,0,lty=3,lwd=1)

axis(1,at = seq(0,1,.2), label = seq(0,1,.2), mgp = c(.8,2,1), cex.axis=1.75)
axis(2,at = seq(-.2,.2,.1), label = rep(" ",times=5), mgp = c(.8,2,1), cex.axis=1.75,las=2)

mtext("75th precentile 1,177 Cases",side=3,line=1.25,cex=1.25,outer=F)
mtext("(e.g., Erie County, New York on April 8)",side=3,line=-0.3,cex=1,outer=F)

mtext("Partisan Animosity",side=1,line=4,cex=1.25,outer=T)
mtext("Marginal Effect of ",side=2,line=5,cex=1.25,outer=T)
mtext("Republican on Policy",side=2,line=3,cex=1.25,outer=T)
dev.off() 

##Figure 7 
pdf(file="figure7.pdf",
    height = 8.5,
    width = 11)
 worrypol.me=c(-0.00256,0.00541,-0.13085,0.00780)
 worrypol.se=c(0.03634,0.03553,0.05031,0.05928)

 behavpol.me=c(1.43914,0.96656,-0.11125,-1.04182)
 behavpol.se=c(0.56482,0.78191,0.69373,0.72292)

policypol.me=c(-0.01570,-0.00096,-0.15898,0.07023)
policypol.se=c(0.04153,0.03730,0.07106,0.08583)

ruler=c(1,2,3,4)

par(mfrow=c(1,3),oma=c(6,3.5,4,1),mar=c(1.5,3.5,1.5,3.5))
plot(ruler,worrypol.me,ylim=c(-.3,.3),xlim=c(1,4),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,worrypol.me-1.96*worrypol.se, 
         ruler,worrypol.me+1.96*worrypol.se, lwd =  1)
segments(ruler,worrypol.me-1.64*worrypol.se, 
         ruler,worrypol.me+1.64*worrypol.se, lwd =  3)
points(ruler, worrypol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,4,by=1), label = c("53 Cases","1,177 Cases","53 Cases","1,177 Cases"), las = 1, tick = T, cex.axis =1.15, line=1.7, col="white") 
axis(1, at = seq(1,4,by=1), label = c("Democrat","Democrat","Republican","Republican"), las = 1, tick = T, cex.axis =1.15)
axis(2, at = seq(-.3,.3,by=.1), label = c(-.3,-.2,-.1,0,1,2,3), las = 1, tick = T, cex.axis =1.7)
mtext("Worry about COVID",line=-1,side=3,cex=1.25)


plot(ruler,behavpol.me,ylim=c(-3,3),xlim=c(1,4),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,behavpol.me-1.96*behavpol.se, 
         ruler,behavpol.me+1.96*behavpol.se, lwd =  1)
segments(ruler,behavpol.me-1.64*behavpol.se, 
         ruler,behavpol.me+1.64*behavpol.se, lwd =  3)
points(ruler, behavpol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,4,by=1), label = c("53 Cases","1,177 Cases","53 Cases","1,177 Cases"), las = 1, tick = T, cex.axis =1.15, line=1.7, col="white") 
axis(1, at = seq(1,4,by=1), label = c("Democrat","Democrat","Republican","Republican"), las = 1, tick = T, cex.axis =1.15)

axis(2, at = seq(-3,3,by=1), label = seq(-3,3,by=1), las = 1, tick = T, cex.axis =1.7)
mtext("Changes in Behavior",line=-1,side=3,cex=1.25)

plot(ruler,policypol.me,ylim=c(-.4,.4),xlim=c(1,4),ylab=" ",xlab=" ",axes=F,type="n")
segments(ruler,policypol.me-1.96*policypol.se, 
         ruler,policypol.me+1.96*policypol.se, lwd =  1)
segments(ruler,policypol.me-1.64*policypol.se, 
         ruler,policypol.me+1.64*policypol.se, lwd =  3)
points(ruler, policypol.me,pch=21,cex=1.5,bg="black")
abline(0,0,lwd=1,lty=3,col="grey30")
axis(1, at = seq(1,4,by=1), label = c("53 Cases","1,177 Cases","53 Cases","1,177 Cases"), las = 1, tick = T, cex.axis =1.15, line=1.7, col="white") 
axis(1, at = seq(1,4,by=1), label = c("Democrat","Democrat","Republican","Republican"), las = 1, tick = T, cex.axis =1.15)
axis(2, at = seq(-.4,.4,by=.2), label = seq(-.4,.4,by=.2), las = 1, tick = T, cex.axis =1.7)

mtext("Marginal Effect of Partisan Animosity",line=1,side=2,cex=1.25, outer=TRUE)
mtext("Support for COVID Policies",line=-1,side=3,cex=1.25)
dev.off() 



